#!/usr/bin/perl -w

###############################################################################
# Transform the geno input files from the platform affymetrix
# into a database load file for the geno table in the database
#
# Parameters needed: dbname, username, geno_files, study, temp_dir
# Output files: temp_dir/geno.out
#
# Author: Herodotos Herodotou and Nedyalko Borisov
# Date: March 11, 2010
#
##############################################################################

use getOpts qw(get_options println);
sub create_ref_mappings;
sub flush_output_cache;

if ($#ARGV != 0)
{
   println qq(Please specify the .properties file);
   println qq(Usage: perl $0 load.properties);
   exit(-1);
}

# Get input parameters
%opt = get_options($ARGV[0]);
$DB_NAME  = $opt{'dbname'};
$USERNAME = $opt{'username'};
$GENO_FILES = $opt{'geno_files'};
$STUDY_ID  = $opt{'study'};
$TEMP_DIR = $opt{'temp_dir'};


# Create the reference mapping for idlinks
$query = qq(SELECT exp_id || ',' || idlink_id FROM idlink WHERE study_id = '$STUDY_ID');
%idlinks = create_ref_mappings($query);

# Create the reference mapping for manufids
$query = qq(SELECT manuf_id || ',' || anno_id FROM anno WHERE platform = 'a');
%manufids = create_ref_mappings($query);

# Clear out any previous output files
!system qq(rm -f $TEMP_DIR/geno_*.out) or die $!;


# Perform the transformation
@genofiles = split(/;/, $GENO_FILES);

foreach (@genofiles)
{
    $genofile = $_;
    $genofile =~ s/\s+$//;
    $num_records = 0;

    # Create a map to cache the output per file (idlink)
    %idlink_to_output = ();

    # Determine the separator based on the file extension
    $data_sep = ($genofile =~ /.*\.csv$/) ? "," : "\t";

    # Open input geno file
    println qq(Transform geno file: $genofile);
    open INFILE, "<", $genofile or die $!;

    # Process the line with the exp ids
    @expids = ();
    while($line = <INFILE>)
    {
       if ($line =~ /^probeset_id/)
       {
          $line =~ s/\s+$//;
          @expids = split($data_sep, $line);
          last;
       }
    }

    # Process the geno data
    $num_line = 0;
    while($line = <INFILE>)
    {
       # Get the data
       $line =~ s/\s+$//;
       @data = split($data_sep, $line);
       $num_data = @data;

       for ($i = 1; $i < $num_data; $i++)
       {
          # Check if exp_id and rsid exist in the database
          if (exists $idlinks{$expids[$i]} && exists $manufids{$data[0]})
          {
             $snpval = ($data[$i] == -1) ? -9 : $data[$i];
             
             # Cache the output
             push(@{$idlink_to_output{$idlinks{$expids[$i]}}}, qq($manufids{$data[0]},$snpval));

             ++$num_records;
          }
          elsif (exists $idlinks{$expids[$i]})
          {
             println qq(ERROR: Unable to find manufacture id $data[0] in the database);
          }
          else
          {
             println qq(ERROR: Unable to find exp_id $expids[$i] in the database);
          }
       }
       
       $num_line ++;
       if ($num_line == 10000) {
          # Flush the output to the files and clear out the cache
          flush_output_cache();
          %idlink_to_output = ();
          $num_line = 0;          
       }
    }

    if ($num_line != 0) {
       # Flush any remaining output to the files
       flush_output_cache();
    }

    # Close file
    close INFILE;
    println qq(Done geno file transformation - produced $num_records records);
}



###############################################################################
# Function: create_ref_mappings
#
# Input: query
#
# Description: This function creates a returns a mapping key->value.
#   When the input "query" get executed, a file in the TEMP directory gets 
#   created. The file must contain lines of the form "key,value". Then, the
#   file is parsed, and the mapping is created.
#
##############################################################################
sub create_ref_mappings()
{
   my $query = $_[0];
   
   # Create the reference file
   !system qq(psql -U $USERNAME -d $DB_NAME -c "$query" -t -o $TEMP_DIR/temp_ref) or die $!;

   # Load the reference file in the map
   open INFILE, "<", "$TEMP_DIR/temp_ref" or die $!;

   my %mappings = ();
   while (my $line = <INFILE>)
   {
      $line =~ s/^\s+//;
      $line =~ s/\s+$//;
      if ($line ne "")
      {
         my @data = split(/,/, $line);
         $mappings{$data[0]} = $data[1];
      }
   }

   close INFILE;

   # Delete the file
   !system qq(rm -f $TEMP_DIR/temp_ref) or die $!;
   return %mappings;
}


###############################################################################
# Function: flush_output_cache
#
# Input: none
#
# Description: This function flushes the output to the geno_<idlink>.out
#              files. Note that it uses the map "idlink_to_output"
#              defined in main.
#
##############################################################################
sub flush_output_cache()
{
   while ((my $idlink_id, my $out_pairs) = each(%idlink_to_output)){
      # Open the file and print the output
      open OUTFILE, ">>", $TEMP_DIR."/geno_".$idlink_id.".out" or die $!;

      foreach my $out (@{$out_pairs}) {
         print OUTFILE qq($out\n);
      }
      
      close OUTFILE;
   }
   
}

